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ABSTRACT 


According to the linearized solutions for thermal 
blooming, the density perturbations become infinite (i.e. 
"catastrophic" defocusing) as the Mach number approaches 
aniey: However, the nonlinearities in the transonic equations 
cutoff the trend to infinity, and the values of the flow 
perturbation quantities are finite. The nonlinear equations 
with heat addition are transformed take Saapis linear 
algebraic equations through the specification of the 
streamline geometry in the heat release region. At a Mach 
| number of unity, streamtube area variation was found to be 
directly proportional to the change in total temperature. 

A steady, two-dimensional mixed flow solution has been found. 
for the transgenic thermal blooming problem. The solution 
for the density perturbations within a laser beam at a Mach 
number of precisely unity is given. For a Gaussian bean 
with an intensity of 3.333x10/ Watts/m* and an atmospheric 


absorption of 8.0x107/ om the maximum fractional density 


perturbation is 1.028x107°, The transonic thermal blooming 


problem does not pose as serious a problem as previously 


anticipated. 
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NOTATION 


area, square centimeters 

speed of sound, meters/second 

function, dimensionless 

line contour around integration region R 
speed of light, meters/second 


specific heat at constant pressure, 
Joules/kilogram mass - degree Kelvin 


specific heat at constant volume, 
Joules/kilogram mass - degree Kelvin 


function, dimensionless 


Green's function, dimensionless; also specified 
function, dimensionless 


Heat quantity, dimensionless 
specific enthalpy, Joules/kilogram mass 
heat addition function, Watts/cubic centimeters 


laser beam intensity, Watts/square centimeter; 
also integral function representation 


unit function; zero for x < 0 and unity for x > 0 
Gladstone-Dale constant, cubic meters/kilogram mass 
scaling factor, dimensionless; length, meters; 

and functional representation of normalized velocity, 


dimensionless 


characteristic length (subscript 1 for x-direction 
and 2 for y-direction), centimeters 


Mach number, dimensionless 


Number of waves, dimensionless 


normal coordinate, centimeters; also index of 
refraction, dimensionless 





Pr 


Re 


Greek 


Static pressure, Newtons/square meter 
Prandtl number, dimensionless 
heat source, Watts/cubic centimeter 


heat source, Watts/square centimeter; also, 
strength of line heat source, Watts/centimeter 


universal gas constant, Joules/kilogram mass - 
degrees Kelvin ; radius of curvature, centimeters; 
and region of integration 

distance, centimeters 


Reynolds number, dimensionless 


line contour around shock waves; also, distance 
along a characteristic, centimeters 


specific entropy, Joules/kilogram mass - degree 
Kelvin ; also streamline coordinate 


static temperature, degrees Kelvin 
maximum streamtube thickness, centimeters 
velocity vector, meters/second 

speed (x-direction), meters/second 

speed (y-direction), meters/second 


laser beam spot size, centimeters; 
also speed, meters/second 


coordinate, centimeters 
coordinate, centimeters 


streamtube shape, centimeters 


atmospheric absorption, l1/centimeters 


variable, dimensionless; (1-4) 172 


flow, dimensionless; 
and (M* -1)2/2 


for subsonic 


for supersonic flow, dimensionless 


ratio of specific heats, dimensionless 


10 











small change, dimensionless 


6 slope, dimensionless; also delta (impulse) function 

€ dimensionless small quantity, 0 < e« << l 

n absolute viscosity, kilogram mass/meter-second; 
also coordinate (y-direction), centimeters 

8 flow angle, beeen 

K thermal conductivity, Watts/centimeter - 
degree Kelvin 

A wavelength, centimeters 

u Mach angle defined by arcsin(1/M) 

v inverse radius of curvature, 1/centimeter 

E coordinate (x-direction), centimeters; also 
transonic similarity parameter, dimensionless 

fe) density, kilogram mass/cubic meter 

z summation 

fo} standard deviation of Gaussian beam, centimeters 

T thickness ratio, dimensionless; also, vibrational 
relaxation time, seconds 

$ viscous dissipation function; Watts/cubic centimeter 

0) velocity potential, dimensionless 

w viscous force vector, Newtons/cubic meter 

Q angular velocity, radians/second 

Superscripts 


( )' 


() 


(~) 


prime, dimensionless quantity 


bar, normalized quantity; also average value of a 
quantity 


tilde, perturbation quantity 


al 





Subscripts 


a,b indicate values of a quantity in front of and 
behind a shock wave 


a x-direction or streamline index - first subscript - 
(i =1 denotes freestream conditions) 


5 y-direction or normal index - second subscript - 
(j =1 denotes centerline conditions) 


L linear value 
0 - total or stagnation condition; also, peak value of 
an intensity distribution 
W wall or surface value 
eo freestream condition 
. Ly ele first, second, etc. perturbation quantities 
g Miscellaneous 
: Vv Del operator, 1/centimeters 
v2 Laplace operator, 1/square centimeters 
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I. INTRODUCTION 


The optical quality of a propagating laser beam is de- 
graded when there are density gradients in the medium. A 
laser beam propagating through an absorbing medium creates 
local ehanges in the density which affect the index of refrac- 
tion of the medium. Therefore, the density gradients will 
refract the light into regions of higher index of refraction 
{higher density) causing defocusing of the beam. This pro- 
cess, which is known as thermal blooming, has been studied 
by many authors [Refs. l1, 2, and 3] to determine its effect 
on a laser beam propagating through various media. 

When a laser beam is slued, a complicated coupled pro- 
cess of thermal expansion 294 fluid flow (relative winds) 
takes place. Depending on the rate of sluing, the beam will 
experience subsonic, transonic, and supersonic winds at vari- 
ous radial distances from the laser source. The subsonic 
and supersonic thermal blooming problems have been > 
solved [Refs. 4-8]. A computer program (DIDER) has been 
developed and experimentally verified [Ref. 9] for beam 
interaction uiehreeuereeaie flows internal to a Gas Dynamic 
Laser. 

A solution for the transonic regime, and, in particular, 
the case when the freestream Mach number is precisely unity, 
has been of considerable interest Hecause the addition of 


heat can lead to extremely strong density gradients. As 
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with the theoretical investigation of transonic flows over 
airfoils and wedges, the linear small perturbation flow 
equations valid for subsonic flow (elliptic equations) and 
‘supersonic flow (hyperbolic equations) are no longer complete 
in the transonic regime. 

Many authors [Refs. 10 - 15] have studied various methods 
of solving the transonic flow equations without heat addition. 
The assumptions based on known transonic behavior about air- 
foils and wedges are used to simplify the governing equations 
and, further, to obtain approximate numerical solutions... 
Often, the assumptions greatly restrict the types of transonic 
flows that can be solved. For example, in the parabolic 
method [Ref. 16], the transonic flow equation is transformed 
into a parabolic differential equation. Since the solutions 
to parabolic differential equations are well established in 
-the literature, this method may be used when the acceleration 
or deceleration is approximately constant in the region of 
interest (i.e. across the cord of many airfoils); however, 
this method is not valid in the region of the stagnation 
point or where the flow acceleration changes sign. 

With heat addition, the basic flow equations for transonic 
flow become more complicated. Zierep [Ref. 17] has solved 
the case when the Mach number is precisely unity; but, in 
solving the equation, he has reduced the problem to that of 
one dimension. Others [Refs. 18 and 19] have solved the 
nonsteady transonic equation with heat addition; however, 


the results are valid only for a short time after beam turn 
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on (i.e. no steady state solution). Even for the short time 
solutions, extrapolation techniques between the property 
values at the critical points on the subsonic side and the 
supersonic side are employed to obtain the values across 

the transonic regime. Experimental investigations are being 
conducted [Ref. 20] to determine the influence of transonic 
sluing on thermal blooming. These experiments indicate that 
a steady two-dimensional solution is possible. Further 
experimentation will provide better guidelines for a 
theoretical analysis of the problem. 

In this thesis, transonic thermal blooming in steady 
two-dimensional invisid flow is formulated and solved for a 
Gaussian heat release distribution. In particular, it examines 
the problem when the freestream Mach number is precisely 
unity. This is accomplished by using a method developed by 
Broadbent [Refs. 21 - 25] in which the streamlines through- 
out the flow are adjusted in such a manner that the heat 
addition necessary to achieve these streamlines equals the 
required heat distribution and boundary conditions. The 
method is exact since specifying the streamlines reduces 
the nonlinear coupled momentum equations in natural (stream- 
line) coordinates to simple algebraic linear difference equa- 
tions that can be "marched" through the entire flow field. 
Boundary conditions are the freestream properties upstream 
of the heat release region and those on a bounding stream- | 
tube sufficiently removed from the heat release region that 


the flow is considered isentropic. The boundary conditions 
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on the bounding streamtube are calculated by methods similar 
to those employed in determining the velocity and pressure 
distributions over airfoils of known shape. Since stagna- 
tion points do not exist on the bounding streamtube flow, 

a modification to the method presented by Sprieter and 
Alksne [Ref. 10] for airfoils and wedges was developed. 

In Section II, the theory used in analyzing the transonic 
thermal blooming problem is presented and includes a des- 
cription of the assumptions made during the mathematical 
analysis, a summary of the pertinent equations derived in 
the Appendices, and a discussion of the method of solution. 
In Section III, the numerical results for the Mach 1.0 flow 
with a Gaussian heat release distribution are presented. 
Lastly, Section IV gives a brief digeusaion of the conclu- 
sions and of the other methods that show promise in solving 


the transonic thermal blooming problem. 
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II. THEORY 


A laser beam that is slued through the atmosphere at a 

. constant rate of 2 radians per second will experience regions 
of subsonic, transonic, and supersonic flow. in planes per- 
pendicular to its axis at various radial positions along the 
beam, Fig. 1. The region to be investigated in this paper 
is the transonic region in which mixed subsonic and super- 
sonic flow exists, as in contrast to purely subsonic flow or 
‘@uversonig flow, and, in particular, the case when the Mach 
number is precisely unity; that is, at a radial position 


given by Eq. (1). 





rE = r*¥ = 79) (1) 


For this analysis, the laser beam is initially of a given 
axial symmetric Gaussian intensity sluing into an isotropic, 
quiescent medium of constant aevsity op (kgm m™ >). Viscosity, 
and thermal conduction are neglected. Justification for 
this assumption is given in Appendix A. The heating effect 
of the laser beam on the medium is approximated by a molecular 


relaxation process which is rapid enough that the absorbed 





energy is instantaneously transformed into heat. The energy 


distribution is of the form given in Eq. (2). 


eee | 
a. exp(- ace ae B : 7) ) (2) 
O 


Q=Ia= 1, 
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where Q (Watts om >) is the time rate of heat absorbed per 
volume; I, (Watts m*) is the beam maximum intensity; 
1 (em 4) is the absorption coefficient; and the exponential 
of the Gaussian with standard deviation o. (cm), which is 
related to the spot size of the beam by w = fo, is 
centered at the origin of the x-y coordinate system trans- 
verse to the beam axis at re, | 

The theory employed for this investigation separates 
the steady transonic flow with heat addition into two en- 
tirely different flow regimes, Fig. 2. The internal, heat 
addition region can be treated exactly by the Broadbent 
method [Ref. 18], in which the nonlinear equations of motion 
and the energy equation are solved numerically in natural 
coordinates. Then, a mesh of streamlines and normals is 
constructed throughout the flow field, Fig. 3., which reduces 


the momentum equations to simple algebraic equations solvable 





for the flow properties. Finally, the mesh of streamlines 
is adjusted until the properties calculated from the linear- 
ized momentum equations satisfy the energy equation and 
boundary conditions. Boundary conditions. are obtained from 
the external, isentropic flow region by approximating the 
transonic solution over a bounding streamtube to the internal 
' region. 

The governing equations in the internal region have been 
derived for flow with heat addition, but without viscosity, 


heat conduction or body forces. A small perturbation analysis 
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for transonic flow (Appendix A) shows that viscosity, heat 
conduction, and body forces can be neglected to first order. 
The small perturbation, steady two-dimensional transonic 


equation with head addition is written as: 


2 x: = tt ag(x',y') 
(1-M Oty Foy ty WM (y+1) Py Pyryt = Cpt, xx 
(3) 


where the term on the right-hand side is added for flows 
with heat addition; q(x',y') (Joules kgm) is the heating 


1 ox-1) is the heat capacity; 


function; Cp (Joules kgm 
T. (°K) is the freestream static temperature; Mis the 
freestream Mach number; and » is the nondimensional velocity 
potential. Therefore, in vector notation, the steady state 


equations [Ref. 26] are: 


V°(pU) = 0 (4) 
pU-vU + vp = 0 (5) 
pU'v(h + U'U) = Q = Ia (6) 
p=oRT | | | (7) 


where U (m sec +) is the velocity; p (newtons m2) is the 
static pressure; h (Joules kgm +) is the enthalpy; and 
R (Joules kgm > °K) is the universal gas constant. Using 
a two-dimensional natural coordinate system [Ref. 26], 


shown in Fig. 3, the governing equations were derived from 
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the above equations. The momentum equations are represented 


in the following convenient numerical form (Appendix B). 





AN (4 441) »3 @ a ae ee 
M Pia, jt+1 i,jtl 
ere (8a) 


~ 9 ee ee 
Yi¢2,g+1 7 “L,54+17 


ym? 


x Lt si ' ' ' ' a 
Bien gel “PL, jeu * oben, 542¥ien, ger 7 Ue, 4% der, 5) 79 


(8b) 


where Ani is the streamtube width at the i,jth station; 

y is the ratio of specific heats; Yi,4 is the inverse curva- 

ture defined as Any 78 and the prime (') denotes dimension- 

less quantities wae respect to the freestream values. The 

3 quantities with a bar and parenthesized subscripts represent 
average values between the indicated stations (i.e. ANG i+) ,j 

is the average streamtube width between station i,j and i+l,j), 

and the properties with the tilde (~) denote perturbation 

values from freestream conditions (i.e. ts = (uss = uy, 5)? ; 

The natural coordinates s and n are tangent and normal to 

the streamlines, respectively. Eqs. (8a) and (8b) represent 

two nonlinear equations since t', §', An' and v' are dependent 


variables. However, specifying the streamlines of the flow, 


these equations become linear equations solvable for 





and 6: at point e in terms of it! 


=f 
if1;g” “itt,3 
at point b, Fig. 3. 


ti! 
itl,j+l i+1,j+1 


at point d and re 441! Bi, j+1 


The remaining flow properties 6: and T! can 


iad j+1 itl,jtl 


be calculated from the continuity equation, Eq. (9), and 
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. equation of state, Eq. (10), derived in Appendix B. How- 
ever, the new properties at the point e (i+1,j+1l) may not 


Satisfy the energy equation, Eq. (11). If it is not satis- 





fied, an iteration process is performed by varying the 
streamline distribution and recalculating the properties 
until coincidence exists throughout the entire flow field. 
Hence, the solution of the problem is reduced to finding 
the streamline distribution that satisfies Eqs. (8) through 
(12) with the boundary conditions determined from the 
external region. 
' 3! ' = ) 
Pi41, jt itl, j+17i+1, 541 ~ + oe 
Pisa, +1 7 Phen, g41 8 he, 542 (10) 


~~ ~ 


Jee. 1, 2 52 ee eae - 
Tita, gta ~ TE, gta 7 Me OY) tied jen 7 UL je2 








a (11) 
(IagAn') (|. ; : 
Oo. (i,i+t1) ,j Asi itl 
where 
2 ,.¥,,CPT,, (1+5-(y-1) M,*) PU, (145 (y-1) M7) ¥ ee 
ag anemccreien : = 
An, 5 AY TET OR fg 


Insight into the behavior of the streamlines at the 
sonic speed was obtained from the solution of the linearized 
small perturbation equation with heat addition [Refs. 4 - 


8]. A computer program (BLOOM) was devised (Appendix G) 
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which gave solutions for a flow very slightly Helou thé sonic 
speed (M. = 0.999), Figs. 4 and 6, and very slightly above 
the sonic speed (M = 1.001), Figs. 5 and 7. ‘Computation 

of the nonlinear term in Eq. (3) from the linear results 

i datcate that the linear equations are valid to M, = 0.9999 
subsonically and M, = 1.0001 supersonically. From the linear 
solutions, a trend of the streamline shape is obtained in 

the near sonic regime, even if not at sonic conditions. 
Nature usually exhibits a strong propensity to accomplish 
changes in a rather smooth fashion (even a shock wave is a 
smooth change in conditions if viewed microscopically 

enough) ; thus, a knowledge of flow conditions in the near 
sonic range gives a good clue to the sonic flow behavior. 

In Figs. 4 through 7, the Gaussian heat distribution is 
centered at the origin on the centerline, and the relative 
airflow is from left to right. 

From the appropriate density calculations and plots, 
discussed later in the thesis, it can be seen that the wake 
term becomes less significant compared with the source term, 
which is directly associated with the horizontal velocity 
perturbation, as sonic jelceiey is approached from either 
direction (i.e. M< 1 and M> 1). Figs. 4 and 5 show the 
horizontal (longitudinal) velocity perturbation at Mach 
numbers of 0.999 and 1.001, respectively. Figs. 6 and 7 
show the vertical (transverse) velocity perturbation at these 


same Mach numbers. Note that, subsonically, the horizontal 
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perturbation, Fig. 4, decreases to a negative valley 
slightly upstream of the laser beam center and then changes 
rapidly to an equal positive peak slightly farther downstream 
of the beam center. The perturbation contours are of very 
small curvature near the beam center; the supersonic per- 
turbations, Fig. 5, show a negative valley near the beam 
center with a return to freestream conditions “downstream of 
the center. In the supersonic case, the contours are nearly 
normal to the flow throughout, in keeping with the character 
of supersonic flow; whereas, in the subsonic case, the con- 
tours are closed, in keeping with subsonic flow behavior. 
The patterns of disturbance, however, appear to be consistent 
as the sonic speed is crossed. 

The vertical or transverse velocity perturbation, Figs. 
6 and 7, show similar behavior in the two speed regimes. 
As expected, the subsonic flow perturbation profiles are 
closed while those for supersonic flow are open to infinity. 
In both cases, the perturbations exhibit an odd symmetry 
about the main flow centerline. For the subsonic case, the 
perturbations are greater than zero to the left of the rela- 
tive flow and less than zero to the right. This is a direct 
indication of the spreading of the streamlines as the heat 
release region is approached by the air flow. Downstream 
of the laser beam, the: streamlines again become parallel 
with the undisturbed flow direction. The supersonic case 
shows the same patterns except that the spreading effect on 


the streamlines is felt infinitely far in the transverse 
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¥ | direction, roughly along the characteristic lines from the 
disturbances. 


The effect of spreading of the streamlines as the beam 





is passed, which can be constructed from these plots of 
longitudinal and transverse velocity perturbations, is thus 
seen to be quite similar immediately below and immediately 
above the sonic speed. Fig. 2 is a sketch of the stream 
lines. It is ekpected therefore that the streamline behavior 
at the sonic speed will not vary radically from this pattern. 
Based on the streamline distribution in Fig. 2, the 
fornewind relationship between area change istveancane width 
variation) and change in total temperature (heat input) was 


chosen, Eq. (13). 


aT 
aS ee ae | 
nom oe ee (13) 


where T is the local total temperature of the flow and 8 
is a variable. One dimensional influence coefficients with. 
area change and heat addition [Ref. 27], Eq. (14), were 


used to determine the required 8 at precisely Mach 1.0. 


am? _ G(x) aa 


‘ (1-M) 


where 


din(T.) 


din (n) ) 
ax 


2 
a= + (1+yM") 


G(x) = M?(1#5(y-1) M*) (-2 
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In order for Eq. (14) to remain finite and to be defined 
at Mach 1.0, G(x) must tend to zero in the limit as the Mach 
number approaches unity. If Eq. (13) is substituted into 
Eq. (14) and the limit taken, it is found that it must be 
zero at precisely Mach 1.0 (Appendix C). Hence, the stream- 


tube area change is related to the heat input by Eq. (15). 


an! , ae 
med = F(yt1) eed (15) 
1,j ee ; 


where the maximum difference between ni ; and ny 5" and T' 
i ee i,j 
6 ti 


and T! is of the order of 10 ° and, hence, the freestream 


1,3 
values can be used in the denominator. Using this approxi- 
mation for the streamtubes, an explicit formulation of the 
streamtube shape, at any position Si | along a streamtube, 
s 


was obtained (Appendix C). 


1 me 2 . 
t = Emad ] e 
Nia x(Ane 4th; ae came (16a) 
ee De ey 
ta ea ge age ee ee 


where for a Gaussian heat input at Mach 1.0 (Appendix C), 





(y+1) f2o'I a sie s! 
An; _ = 5 exp (-U=F)-) lterf(——)} ; i cae Ses eae 
J . Ss 20! [20! 
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Once the streamtube shape is known, the slope (dy/dx), 


Eq. (17), and the radius of curvature (R), Eq. (18), can 


“w be calculated at every mesh point throughout the flow. 
A 
/ 
= ' ae -n! i ' 
ay = : 44, 4 M5414 nig t 7875. dott 
Ax! ; . .. i,j AS; . . ‘ —u! 
Le] oJ | j>1 
(17) 
v! . = aa ee = ie i-l,j x 93,3 ~ §4-1,5 eee eee 1 
te) Re pet eee aye es ASi 4 i 
J i,j ’ J >1 


Therefore, the proportionality between the streamtube 
i shape and heat release distribution gives not only an 
explicit formulation for the streamtubes and, hence, linear- 
. ization of the nonlinear momentum equations but also simul- 
taneous agreement between those properties calculated from 
the momentum equations and the energy equation. The problem, 
in the case of precisely M. = 1.0 flow, becomes one of 
| determining the correct boundary conditions. 
The boundary senaieiona that are used are those corres~ 
B15 Tl =0), 


which are far enough upstream that the properties are not 


ponding to freestream conditions (ij . = pi . 
“Ld 1,3 


affected by the perturbation caused by the heat release 
region, and the pressure (Pi as) and/or velocity (up se) on 
the j*th streamtube which is considered far enough removed 
from the heat release region that the flow can be considered 


isentropic. Eqs. (16b) and (16c) can be used to determine 
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the streamtube shape and displacement where j = j*. However, 
a more convenient form is derived in Appendix C for a 
Gaussian heat release distribution at M, = 1.0 which can be 


written as: 








(y+) Toao 'n _ 8 . . 
ae aie AGS te ae) eee 
(19) 
and the slope of the bounding streamtube: 
an! .,  (ytl)Ia fare! s!° ieaee 
agree = ——ag ew - eg = ep - (20) 
ind 20! 2710 ' 2c' 


where t' is the maximum growth in the streamtube. 

Since the mass flow is constant ina gLegeteune: the 
bounding streamtube can be considered as a solid wall simi- 
lar to the surface of a wedge with a cusp extending forward 
to the freestream conditions upstream at infinity, Fig. 2. 
By hypothesis, this flow is isentropic and can be treated 
by methods similar to those used in analyzing transonic 
airfoil and wedge flows. In particular, an approximate 
solution is obtained through an iteration process on the 
integral equation derived by Green's functions for small 
disturbance transonic flow similar to that performed by 
Sprieter and Alksne [Ref. 10 andl1d] for thin airfoils and 
wedges but modified for the streamtube configuration which 


has no stagnation points. The absence of stagnation points 
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: in the case of the flow around the bounding streamtube 
enhances the accuracy of the approximate numerical tech- 
niques because the stagnation point represents a singularity 
in the problem and places added constraints on the method and 
accuracy of the solution in the region of the stagnation | 
point. In this method, the quadratic, nonlinear nature of 
the governing equations is retained which precludes shock- 
free supercritical flows in the transonic eee 

The integral equation for transonic flow is derived by 
a Green's function analysis on the small perturbation 
transonic equation without heat addition, Eq. (3) without 


the right-hand side, resulting in Eq. (21). 


— (> ast m2 = py _ 
u(x, 0) = ae + ZU, (x, 0) si (21) 
where 

= mM *(y+1)G'(x,0) mM? - m ? 
UO) (22) 
w B2 37 
—_- = Bee ~ az dé 

(x,0) == ‘oa — (23) 
“Dy Pk ee (x-E) 

i ud (X, 0) F-X es 
r= s 4 — e(22*) af, (24) 
os gE (&7-*) ee Be . 9° 1 a a= 24 
D ag gee eg (24) 





° and the normalized variables are defined in Appendix D. 
Eq. (21) is a function only of x, the normalized coordinate 
along the wall surface, and Tt, the normalized amplitude of 


the wall slope, given by: 





= a, | 
c= 7 exp - (25a) 
x | o | 





where | 
~ Mo? (y+) 1 M*(y?=1) 10 Jor Oo 
T= 3 oa (25b) 
a B B° pi u.yY 


The iteration technique (Appendix D) is performed for 
° various subsonic freestream Mach numbers starting at the 
initial eupereeteical flow and proceeding toward the Mach 
1.0 flow. Each solution obtained by the iteration technique 
gives a unique t, which determines the freestream Mach number 
and its corresponding normalized velocity distribution 


4/3 is plotted 


U,, (X, 0) along the surface. If (a -1) versus T 
for each position along the surface, it is found that as 

the freestream Mach number tends toward unity the slope 
becomes constant. This corresponds to the phenomenon of the 
Mach number freeze [Ref. 26] wherein the local Mach number 


is invariant with changes in the freestream Mach number 


when the latter is near unity or, more precisely, 


= 0 (26) 
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Vincenti and Wagoner [Ref. 28], Liepmann and Bryson [Refs. 
29 and 30], and others have shown that the corresponding 
approximate relation yielded by the small-disturbance 


transonic theory is given by: 


d : | 
ae = Q (27) 
lB 
where 
(1-m?) (G-1) 
== ne = (28) 


(m2 (y#1) 1) 273 


The freeze of Mach number extends over a finite range near 
Mach 1.0 where it is constant. Hence, €, the slope of the 
(u-1) versus T curve at each location x, is constant near 
the sonic eehaeichs aha the local Mach number at each 
position can be calculated by Eq. (28). With the local 
Mach number distribution known, the velocity perturbation 
can be calculated from Eq. (22), and the other flow proper- 
ties can be determined from isentropic flow relations [Ref. 
27). 

The Mach 1.0 flow is completely specified. Boundary 
conditions are known upstream of the heat release region 
and on a bounding streamtube. Calculation of the entire 
flow field is readily calculated from the simple algebraic 
equations obtained by Broadbent's method, where the stream- 


tube configuration is known precisely at Mach 1.0. 
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III. NUMERICAL RESULTS 


The computed flow properties are for the Mach 1.0 thermal | 


a blooming problem with the characteristics listed in Table I. | 









Gaussian Beam 






Peak Intensity (I,) 3.333x10! watts/m” 










5 cm. 


Atmospheric Absorption (a) 8.0 x 107” om? 


(co) 


Standard Deviation 









Atmospheric Freestream 
Conditions 













°K 
5 


288.0 





Temperature (T,) 









Newtons/m* 





Pressure (p,) 1.013x10 










Density (p,) 1.2246 kgm/m> 





Velocity (u, = a*) 340.0 m/sec 







Streamtube Width 
(An, 4=An,,) 
a 











Streamline Cell Length 
(As; 4) 


l om. 
Specific Heat at Constant 
Pressure (Cp) 1005.0 Joules/kgm°K 


Table I 










Freestream Flow Properties And 
. Laser Beam Characteristics 
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The first step in the procedure is to determine the 
boundary conditions for the heat release region. Far up- 
oar | stream of the heat release region, all the dimensionless 

perturbation quantities are zero. The other boundary con- 
ditions are calculated by Sprieter's method [Ref. 10] for 
the bounding streamtube, the shape of which depends upon 
the heat release distribution and is known, Eq. (16). 
Using Sprieter's numerical iteration procedure on the 
| transonic flow integral equation, Eq. (21), the normalized 
velocity distribution along de boundins streamtube “ae 
- obtained for several freestream Mach numbers between the 
initial supercritical flow and Mach 1.0. The results are 
graphed in Fig. 8. From these solutions, a (U-1) versus 
72/3 curve, Fig. 9, was plotted for several positions along 
: "the streamtube. The results indicate that the slope (é&) 
is constant for each position and, hence, the Mach number 
can be considered frozen. | 
From the values of ED), the local Mach number variation 
ata freestream Mach number of unity was calculated from 
Eq. (28) ’ Since — is known and constant for the known 
physical streamtube shape, Eq. (30), from Eq. (D6). 


6 


t= = 1.386 x 10. (30) 


{aro 





where 


| 2 
- (y-1) Lato = “5 
t YD_u,, 1.737 x 10 cm. 
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(a) Shock Wave at x= 9.5 | 
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(b) Shock Wave at x = 11.5 
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“ (c) Shock Wave at x = 14.5 
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(da) Shock Wave at x = 19.5 
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(e) Shock Wave at x = 24.5 
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FIGURE 8. High Subsonic Supercritical Flows on 
- Bounding Streamtube 
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is the maximum growth of the bounding streamtube. 

Once the local Mach number was determined, the velocity 
bel turbaeion along the bounding streamtube was calculated | 
from Eq. (22) and is plotted in Fig. 10. The remaining 
flow properties were calculated by the isentropic flow 
relations [Ref. 27]. | | 

It has been shown that the area change of the stream- 
‘tubes in the heat addition region is directly proportional 
to the change in total temperature; that is, 8 = 0 in 
Eq. (13). The energy equation is automatically satisfied 
by this selection for the area variation. Since the heat 
release distribution is known, all of the streamtubes are 
specified. From the streamtube shape, the radius of curva- 
ture of the streamlines.:were calculated at every point 
throughout the region, Eq. (18). 

With the curvature known, the calculation of the flow 
properties in the heat region is reduced to the solution 
of two simple linear algebraic equations, Eqs. (8a) and 
(8b), the momentum equations in natural coordinates. The 
results of this "marching" technique are shown in Fig. ll 
where the density perturbations in the upper half of the 
symmetric laser beam are presented. This method has the 
distinct advantage that it can be solved on an HP 9830 
computer. 

Figs. 12 and 13 show the linear supersonic and subsonic 
density perturbation solutions calculated from the computer 


program (BLOOM) developed in Appendix G. The results 
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. | presented in Fig. 10 at Mach 1.0 are consistent with a 
continuous transition through the transonic regime. Although 
aes - the solution calculated is for precisely Mach 1.0, it can 
be logically hypothesized that the flow behavior through 
the transonic regime is as shown graphically in Fig. 14. 
The linear siverscuic and subsonic solutions, and 
Sprieter's high subsonic solutions are in excellent agree- 
ment with this hypothesis as well as that sein eetured from 
the hodograph plane representation for transonic flow 
(Appendix E). Therefore, there is a steady, two-dimensional 
ie | solution for the transonic thermal blooming problem which 
is finite and does not result in "catastrophic" defocusing 


of the laser beam. 
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IV. CONCLUSIONS 


It has been demonstrated that there is a two-dimensional 
steady state solution to the transonic thermal blooming 
problem for a sluing laser beam with a freestream Mach number 
of precisely ma, According to the linearized solutions 
for thermal blooming, the density perturbations become 
infinite as a Mach, number of unity is approached... Due to 
nonlinear effects the trend to infinity is cutoff, and the 
finite value of perturbation quantities has been established: 

In obtaining the solution, it was found that at precisely 
Mach 1.0 that the streamtube area variation through the heat 
Seiéane region is directly proportional to the change in 
total temperature in that region. With the streamtube 
shapes known, an exact solution for the flow field can be 
obtained by Broadbent's method [Ref. 21] with the corres- 
ponding boundary values for Mach 1.0 isentropic flow far 
from the heat release region. There are no restrictions on 
the existence of shock waves nor on the symmetry or asymmetry 
of the heat release distribution, provided the bounding 
streamtube in the isentropic region and the boundary 
conditions on it can be calculated. | 

Furthermore, in obtaining the solution at precisely 
Mach 1.0, it was found that shock waves exist throughout 
the transonic regime, Fig. 14. As the sonic condition is 


approached from subsonic speeds, a shock wave forms at the 
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(b) Supersonic Supercritical Flow _—_ 
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(c) Sonic Flow — 


(£) Mach Number (M=l+e) versus Position Along 
Beam Centerline in the Transonic Regime. 
(Dashed line represents hypothesized 
Supersonic behavior.) 


FIGURE 14. Laser Beam Centerline Mach 
Number Distribution in the 
Transonic Regime 





(e) Subsonic Flow 





downstream 1/e (one spot size) size, of the Gaussian heat 
release distribution, where the local Mach number first 
becomes sonic, and rapidly moves downstream to infinity. 
If the sonic speed is approached supersonically, a ahoek 
is formed at the beam center. To be consistent with the 
results in the subsonic portion of the transonic regime, 
it has been hypothesized that the shock moves rapidly up- 
street to infinity. Finally, the analysis in Appendix F 
indicates that there is a minimum laser sluing rate that 
must be maintained to preclude significant phase distortion | 
in the beam at Mach 1.0. 

The results of the Mach 1.0 thermal blooming problem 
for a Gaussian intensity distribution are in excellent 
agreement with the approximate results obtained by Ellinwood 
and Mirels [Ref. 19]. Fig. 15a shows these results where 
the maximum density perturbation was calculated to be of 
the order of 92/3, Q is defined as: 

(y-1) Tar 7 
CF —_—_— _, (31) 
ypa 

where r is the radius (spot size) of a circular beam of 
uniform heating intensity (I). irda appro uaciony when 
applied to the Gaussian beam of Section III, Gives a maximum 


4 which is consistent 


density perturbation of the order of 10° 
with that obtained in this paper, and shown in Fig. 15b. 
Although the density perturbations are very small (i.e. 


measured in parts per million), significant optical defocusing 
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FIGURE 15. Schematic Representation of the Steady State 
Density Perturbations for a Sluing 
Two-Dimensional Laser Beam 
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of the laser beam could result (Appendix F). With an 
increase in the beam area due to defocusing, a large decrease 
in the beam intensity I (Watts/m-) will occur. 

Two other methods that show promise in solving the tran- 
sonic thermal blooming problem are the finite difference 
method developed by Murman and Cole [Ref. 12] for isentropic. 
transonic flow which can be extended to include heat addi- | 
tion, and the "finite element method" of solving partial 
differential equations [Ref. 31] which is a powerful method 
for solving a wide range of engineering problems. 

A Mach 1.0, steady state seiabien for thermal blooming 
ina laser beam with a given Gaussian heat release distribution 
has been presented in this thesis, Fig. 11. Subsequent inves- 
tigations should be performed to ascertain the Mach 1.0 
solutions for various Gaussian and non-Gaussian heat release 
distributions, beam intensities (I), atmospheric absorption 
coefficients (a), and freestream flow properties. In this 
manner, a correlation might be developed that would be use- 
ful as an engineering approximation to the density perturba- 
tions at various transonic Mach numbers, vice a simple order 
of magnitude linear extrapolation as presented by Ellinwood 
and Mirels [Ref. 19]. Although solutions are known at Mach 
1.0 and in the linear subsonic and supersonic regions, from 
the computer program BLOOM, it would be premature to estab- 
lish a criterion that could be used in estimating transonic 


behavior in general before additional work is completed. 
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APPENDIX A 


% : DERIVATION OF STEADY TRANSONIC FLOW 
FOR A CONDUCTING GAS WITH HEAT ADDITION 


Eqs. (Al) through (A4) govern steady compressible flow 
with heat addition and constant properties [Ref. 32]. These 
equations include viscosity and heat conduction. The 
distribution of the heat addition is given by the function 
Qof (x,y) and is assumed known. The gas flow is assumed 


uniform at infinity (x = +). 








V°(pU) = 0 : | (Al) 
; pU-v0 + Vp = y (A2) 
where yp = nV(V-°U)/3 + nv7u 
oTU-Vs = - O- Veg (A3) 
where  o= an (v0) 2 - a (U-U) - (vxt)* - 20-70) 
2 
and Veq =-KV T - Qf (x,y) 
ete 2 
Vp = pa Vs/Cp + a’ Vp (A4a) 
or 
x s-S 
exp(—e—2) = 8/9 V1) (A4b) 
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The following dimensionless variables are defined: 


* 
I 





| = x/hy y' = y/L, 
* o' =' p/p, T' = 7/7, wu = 0/0, 
. s' = s/Cv,, n' = N/ Neo = 1 K' = K/K, = 1 


where Ly and L, are characteristic lengths in the x and y 
directions, respectively, and the reference values are the 
freestream conditions. Substituting the dimensionless 
variables into Eqs. (Al) through (A4) and eliminating the 
pressure term in the momentum equation, Eq. (A2), with the 


equation of state, Eq. (A4), the equations simplify to the 


following: 
= o'V'-U' + U'-Viep' = 0 (AS) 
‘ a. = totnt : ota ' '.q! 12a 
o'u'-V'tt = - an ye Sean Be + a + — (A6) 
M., YM, | 
- y(y-1) M,? 2 Pare: 2 sr) 
o'T'u' V's! = —— { = 5 (Vv! *U') + V! (U'-U) 
- (v'xt')? - 2(8'-9'7a") | 
yv'em 2 
+ Dae tt YD Hf (xy) set} 


where the dimensionless ratios are defined as follows: 


: Reynolds number: Re = 9 UL L,/n., 
Mach number : MS = U,/VyRT 

e Prandtl number: Pr = Cy Neo Kes 
Heat quantity: q.. = 0,2;/0,.U,,” 

me . Scaling factor: L = £,/L5 








For flows in which conduction and eiaeogiey can be 
neglected, two flows will be dynamically similar when y, 
* | M, (Mach number in the absence of heat addition) .and H, 
are the same. | 
Retaining the viscous, aonduction, and heat addition 
terms, the governing equations, Eqs. (A5) through (A7) , in. 


two-dimensions become: 


1,ou! ov! + op! a ae 
esr tig tu xr t BV Gyr 0 (A8) 


x~direction momentum equation: 




















me | | 2 
du' ou! T*- op! p'T' os't 4 du! 
a o'(ut + Lv' =) =- —= os er tl os 
: aX oy mM? ox! YM,” ox 3Re 9x? 

2 2: x2 
L ov' L* ou! 
; SRE OxTby" RE 372 ) 
y-direction momentum equation: 
syegaee: Oe aa _LT' a0' _ bott' ast | 4n? 92yi 
t T “a. Tt ee ¥ 
ox Pe) mM,” oy ym, dy Re ay "2 
. L 9*u! 1 o“v! (A10) 
3Re ox oy’ Re 9x17 
energy equation: 
# (y-1)M 
Meike AOS 1 os’, _ a 4 ou’, 2 _ gu'gv' 
p'T'(u Sxr t Lv ay7 ) Re § (xr ) L axtoyr ay 
a Fs jv! } ei v' +. guiy2 (2 Sa £ L2 P) Ty 
' t 
yy"? ox dy. RePr f2 5y" 
dds 
y(y-1)M,7H,, £ (x,y) | (All) 
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Since the transonic flow region is of interest, a 


lle: 
transonic expansion procedure similar to that used for 
af a incompressible flow around an airfoil of thickness ratio T 
{[Ref. 33] will be used. In the following analysis, the 
heat addition will cause the perturbations from the free-. 
stream conditions. Introducing the following perturbation 
quantities: 
were et are A/a ; vi = ev} 
p'=i1l+t gare p; + 4/3 5 ; T' 1+ 7/3 T, + e473 Ts 
“ 
staree/3 sg: ; a/Re~ Ole) + H,~Ole*/9) | 
‘ = : alfa ts a ta73 
. the first and second order equations at M, = 1.0 become: 
First order equations: 
dus op4 
1 1 ae 
Jur t Pee 0 (Al2) 
8V7 963 
3xr * dy? = 0 (A1L3) 
ats | att 
wert (y-1) ger = 0 (A14) 
og Integrating Eq. (Al2) gives uy = “67 which, when substituted 
e . 
i into Eq. (Al13), illustrates that for first order theory the 
flow is irrotational. 
«x 
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bee ek Seer, | | ~ (A115) 


Finally, Eq. (Al4) gives the temperature distribution as: 


= (y-1) BL = -(y-1) OY (A16) 











ats at; = av] BBS 904 
' : “= 
i + PF yur + byt + *yr + Uy ¥yr = 0. (A17) 
~s ~ “4 ~t ~ 1 
“ au, ae dus a dus we: pa er aa 
ox’ 1 Ox! 1 Ox! M a Ox! 1 0x 
; | a8: 4et/3 92a: 
| YM a 
«. vi vi P5 p 
r avs a ova ae ies 904 ae mee 
1 Ox" * “1 Ox" m2. 2 dy" ay" 
a ee 375s 
-—5 + ¢€ (A119) 
ym, ye ax! 
es 2 
as! oct! | 
—* = yel/3 1, Y(y-1) M : £(x,y) — (A20) 
ox ax! . 


Combining Eqs. (Al17) through (Al19) and using Eqs. (Al2) and | 


(A1l6), the transonic perturbation equation is obtained. 








oul ov! 1/3 o°u r) zs 
(y+1) Uy ax! + y' 3 Qx'2 + 5a" + (y 1) £(x,y) 
> | | (a21) 
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The first term on the right-hand side of Eq. (A21) 
represents the effect of viscosity, the second term the 


‘ | conduction, and the last the heat addition. As can be seen, 


1/3 


both the viscous and conduction terms contain a ¢€ term 


1 


| when the Reynolds number is of the order of € ~, and, hence, 


are small and can be neglected under most circumstances. 
If, however, the Reynolds number is of the order’ of @~ 2/3 

or smaller, both terms become important and must be retained 
in the perturbation equation. Since to first order the flow 


is irrotational, the velocity potentials uy = 0¢/dx' and 


Vi = 0o¢/dy' can be introduced into Eq. (A21) giving: 


1/3 
a (VFL) Or dry Ft oytye FO (3y+l)¢€ _— 


+ (y-1) £(x', 6 t/3y!) (A22) 


A similar analysis can be carried out for a freestream 





Mach number close to unity. In this case, the freestream 


2/3 


Mach number can be approximated by M, = 1 + Ke and the 


transonic expansion gives: 
me (y+1) 91) Oytye — bytyt = 


- (y-1) £(x',e/3y') (23) 


t 
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APPENDIX B 
BRIEF DEVELOPMENT OF BROADBENT'S METHOD OF 
SOLVING NONE SEAS EQUATIONS’. OF MOTION WITH 
HEAT ADDITION 
Broadbent's method [Ref. 24] is an inverse solution in | 
which the streamlines, normals, and boundary conditions 
along at least one streamline and normal are specified and 
the resulting flow field calculated. 


The governing equations for this derivation are as 


follows: 
: V+(p0) = 0 (Bl) 
~ ou-v0 + Vp = 0 (B2) 
. oU-V(h + 5u-0) =Q = Ia : (B3) 
p = pRT | : (B4) 


In natural coordinates [Ref. 26], these equations 


become : 
2_ (puan) = 0 (BS) 
os p 
ou , oP . 

pu ga * oe” °.. (Bo) 
, u , ap | 

ee +55 7 9 | (B7) 
> . 12 
> | pu s- (h + su") = Ta (B8) 
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A mesh is constructed over the flow field using stream- 
lines and normals which are chosen in advance as in Fig. cf 
The cells specify the streamtube boundaries and the mesh 
lines the streamline slope and curvature ehvougneut the 
flow field. With the mesh specified, Eqs. (B5) through (B7) 
< are nondimensionalized using uniform freestream flow 


properties (ay se Py, 4! and Tr 4) a constant uniform 


5 t 
streamtube width upstream of the perturbations (An) 5) and 
a : 


constant Cp. 


t ' 1 st 
Pig a ae 
An‘. . . 
t 4 roe hee i,itl), ~ Ae: 
Bie, j+1 7 UL, 541 age (Pisa, j+1 7 Pi, jen) 
is (B10a) 
p! -p! es, a i YM Zu! yt. 
“itl1,j+l1 it+l1,j 2 co itl, jtl itl,j+1 
+ Ui4en,5%ie1,45) (B10b) 
mt a mt 2 2 2 ant as 
Thea, gen Th, gen TOF Met) ie) 54. 7 b,j 47) 
eee 
es a 1 eae (B11) 
| : i,jtl ~ 21,541 
T! (B12) 


’ ae | 
Pi,j > 94,3 F7i,3 


where Q. is a freestream heat quantity defined in Eq. (12); 


is a dimensionless heating term; and v! . is a 


| 
i,jt+l 1,Jj 


dimensionless inverse radius of curvature. 


64 





The slope (dy/dx) and the radius of curvature (R) of 
the streamlines at a point b (i+l,j) are shown geometrically 
in Figs. 3. Due to the symmetry of the problem, the slope 
and the curvature of the centerline streamline are zero 
(i.e. for all stations i > 1 and j = 1). Assuming that 
AYG = An, and Ax, 3 = As; (i.e. to first order the 
mesh is rectangular), and that the streamline shape is 
known explicitly, the slope in finite difference form is 


expressed as: 





- t tt. 
ay ee Mid, j71/2Ani 15 ie + 1/2 An; 
x i,j Le As; =| 
(B13a) 
or: 
J-1 Lo, Ls | Jt ' eer 
i+, ae {Ah +1,1 FL AL EL 3 ke ij 
1,7) As! : 
ene (B13b) 
where 
nf, = gan! , = 30h! , + 1) (B14a) 
i,t 27 ee 2 tek 
t a jv ’ se jut : 
a cg a ace ae wait G-gr i 21, Bl) 
j3>2 


The inverse of the radius of curvature (vy) is defined in 


Eq. (B15) as: 
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(B15) 





By nondimensionalizing and substituting the slope obtained 
in Eqs. (B13a) or (B13b) into Eq. (B15), the following 


finite difference equation is obtained: 


; | (65) 47 $5.3) 5? | 
An, ;: Kens | 
vi, = —gile in) 5 (B16) 
= (L + (8, 5)°) 


Because the magnitude of the slope S544 is never greater 

than 10°, it can be neglected in the denominator simplifying 
the expression for the curvature. Substituting Eq. (B13b) 
into Eq. (B16) and vataiaine only ectie of first order, the 


curvature becomes 


1 a ae! L.~ g cig ae 
ee [Pita - efia, 4 a ek 
| (B17) 
1 jrl 1,~ 2 
- ant 5 +5Ani 1 ae ani-1,17 sani, | Mas; 4) 


If the mesh is not specified, Eqs. (B9) through (B12) 
represent five equations in five unknowns and are nonlinear 
since An', p' and u' are dependent functions. When the 
streamtube shapes are specified, An' and the curvature y! 
are known quantities and the equations become linear in p' 
and u'. Therefore, Eqs. (Bl0a) and (B10b) represent two 


algebraic finite difference equations solvable for p' and u' 
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eurougheat the £low field, provided suitable boundary condi- 
tions are known. For example, if the values of p' and u' 
are known at points a (i,j), b (it+l,j) andd (i,j+l) in Fig. 
3, the pressure and velocity perturbations can be calculated 
algebraically at point e (it1,j+1). ‘This procedure is con- 
tinued ("marched") throughout the mesh until all flow pro- 
-perties have been determined. | 

The remaining flow properties, density, and temperature 
are then calculated from Eqs. (B9) and (B12), the continuity 
equation and equation of state. Once this has been accom- 
plished, the energy equation, Eq. (Bll), is checked to see 
if the properties calculated from the specified streamline 
shapes agree with the desired heat distribution. If it does 
not coincide, an iteration procedure is per earaed in which the.. 
streamline shapes are varied until agreement is obtained. 

Boundary conditions are needed along at least one stream- 
line (usually a wall) and one normal (usually freestream 
conditions). In the transonic thermal blooming problem, the 
boundary values are chosen as the unperturbed freestream 
properties upstream of the heat release region and the 
velocity or pressure on a bounding streamtube which is far 
enough away from the heat release region that it can ba: 
considered isentropic, Fig. 2. With the pressure and/or 
velocity known on two of the boundaries, the goraeion for 
all the flow properties can be calculated as mentioned pre- 


viously. If a coundary condition is known at other points 
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in the flow field, the resulting properties from the 

: "marching" procedure must agree with it or the specified 
streamtube configuration must be altered until coincidence 
between the boundary conditions and heat distribution is 


obtained throughout the flow field. 
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APPENDIX C° 
DERIVATION FOR THE RELATION BETWEEN AREA 
VARIATION AND HEAT ADDITION AT PRECISELY 
MACH 1.0 FLOW 
1. AREA CHANGE WITH HEAT ADDITION AT THE SONIC POINT 
From the influence coefficients for constant specific 
heat and molecular weight [Ref. 27], the change in Mach 


number with heat addition and area variation is given by: 


2 


aM“ _ G(x) 
= = (Cla) 


where 


ene Lae 


2g ina) 
— > - og), (cb) 


G(x) = MO (145(y-1) M2) [ (L+yM?) 


d (int) /dx and d(lnA)/dx are assumed to be functions of x, 

the streamwise coordinate, and dx is always positive. Based 
on the relation between streamline shape and total temperature 
change obtained from the linear subsonic flow and supersonic 
flow, Figs. 4 through 7, the following relationship between 
area variation and total temperature change (heat addition) 


is chosen: 


d(intv_) 
gine) = (y+) (148 (x) ) —_— : | (C2) 


66 





% 


8 is a variable providing an additional degree of freedom 
in determining the relationship at Mach 1.0. Substituting 
Eq. (C2) into Eq. (Clb), G(x) simplifies to: | 


d(ln T_) 


G(x) = MO (145(y-1)M*) ((M7-1) 7 - (y#1)8) g—S. (C3) 





Since the denominator of Eq. (Cla) tends to zero in 





the limit as the Mach number approaches unity, the numerator 
wee go to zero to permit a continuous passage from subsonic 
to supersonic speeds (or vice versa). The condition, when 
G(x) passes through zero simultaneously with the Mach number 
becoming unity, is designated the critical point G(x) and, 
from Eq. (C3), is equal to: 


2 d(in T,) 


Bas = 0 (C4) 


lim G(x) = G*(x) = ~ 2 (y41) 
M>1 


Hence, for heat addition -d(ln T,) /ax > 0, 8 must equal 


zero or, equivalently, the area change is directly propor- 


tional to the heat addition (i.e. change in total tempera- 


ture) for steady flow from Mach 1.0, Eq. (C5) is: 


: d(ln T_) 
giin A) = (yt) 2. (C5) 


At the critical point, dm? /ax is of the indeterminate 


form 0/0. Therefore, the limiting value of am? /ax at the 


critical point is found by applying L'Hospital's Rule to the 





right-hand side of Eq. (Cla). 


2 Bors 
. GM” _ ,dM (dG/dx) * 
oe oe ~(aM*/ax) * 


M>1 
Or, equivalently, 
(aM?) de? he ~ (2¢ * . (C6b) 
dx ~ ¥ ax a , 


Eq. (C6b) suggests that (am? /ax) * can be positive or negative 





and the flow can become either subsonic or supersonic from 
the sonic point. The requirements on 8 follow from Eq. (C6b) 


where G(x) is given by Eq. (C3). 


| | 2 d(ln T) 
dG Bo poe: HAG? ps2 ake - dB aM re) 
Sag) oe aa (yell ad + Vas 7 7 aaa 

B=0 (C7) 


Thus, at the critical point, the flow may be continuous 
through the sonic speed provided (dG/dx)* is negative and 
is not zero. For heat addition, Eq. (C7) can be solved for 


the critical value of (d8/dx)*, Eq. (C8): 


(C8) 





be > ty 


Y 


and is graphed in Fig. lc. 
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Allowable Solutions for 


Heat Addition d 1n (T,) /dx >0 | “Slope y/(y+1) 
> 
aM", x 
ra 
Supersonic 
Subsonic Allowable Solutions for 


Cooling dn (T,) /dx <0 


FIGURE 1C. 8 Variation versus Mach Number Variation 
for Sonic Flow 


From the sonic flow condition, 8 must be positive for 
supersonic flow and can be either positive or negative 
for subsonic flow. The change in 8 with Mach number squared 
_ must be greater than y/(y+l) for either flow condition, 
supersonic or subsonic. 

In conclusion, for sonic flow the area variation is 
directly proportional to the change in total temperature, 
Eq. (C5), where the change of 8 from the sonic point must 


satisfy Eq. (C8). 
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2. STREAMLINE CONFIGURATION FOR MACH 1.0 FLOW WITH A 
GAUSSIAN HEAT INPUT DISTRIBUTION 


In the first part of this Appendix the streamtube area 
change was shown to be directly proportional to the change 
in total temperature at precisely M, = 1.0; see Eq. (C5). 

To derive an expression for the change in total temperature, 


a slightly different form of the energy equation is used, 


Eq. (C9), 


te] 


pu z= (CpT,) = Q = Ta (C9) 


my 


in terms of the total temperature T,. Assuming a constant 
Cp and a Gaussian heat intensity distribution in the form 


given in Eq. (2), Eq. (C9) becomes: 


a ne 2 
at 
° Ia +2 ss 
7 tyJ = aT f exp - _, dy' exp - 5 ds: 
O4 * co 1,j ' 20' ? 
oJ -n! . 
1s) 
(C10) 

where n; , Tepresents the total displacément of the jth 

v . : 


streamtube from the centerline’ (j=1), and Q,, defined by _ 
Eq. (12), is evaluated at M, = 1.0.. 

Equation (C10) is now substituted into Eq. (C9) to give 
an expression for a change in total streamtube growth over 


a streamline distance As; at MoS 1.0. 





Ont 3 2 Ss: 
+1)1 i 
ee geet, fs Oe We 
2 o l,j » <2" 20' 
eee ) 
tJ (C11) 
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or more conveniently for computational purposes: 

















~ 
t 
. are | 2 
_— (y+1)I_o a iiye Ss: 
' - a O a. ’ _ L nas ' 
| dni 0, | f exp ae dy €xp San ds: 
1 10 s;° 
ot ; . oO = Es t 7 ° — 
q(v+1) . exP so) dsj i>1, j=1 | (C12a) 
and 
o 1 1 
* '.= : + =(Ytl) — 
Sg nigh. 
e j-l (22+1) An, a/2 2 12 
. ’ 1 i 
y exp - i, dy | exp - 5 ds 
af ~(2£-1) (An, o/2) 20! 20° 
v 
1 I j-1 (RAn} 1" at 
= : + =(Y¥+tl) — i, ecciceisl - =e dg! 
dni x (Y 1) o. ) exp eae exp i Si 
L=]1 
(C12b) 
for streamtubes off of the centerline (j 7 1). Equations 
(Cl2a) and (C12b) can be integrated to give an expression 
for the total streamtube perturbation from its initial 
: ; : 
width n1,3 as: 
An! I_av2 7 s! 
' an ee 1iztl fe) i } oa se 
: i,17 "4,1 5 x(¥+1) om {1 tere = a Cae dae B 
| (C13a) 
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ee "= | 
’ eee | = = ' ' 
i,j "143 2 Ani. ' ) ANT 241 

I_av2 oc! 
~ i ~t L ee 
* 92 Si | 
exp - ——z (1 + erf )} 7. /2o1, jal . 
R=1 20! ¥2 o! : 


(C13b) 


An alternate form of Eqs. (Cl3a) and (C13b) which can be 
used in Eqs. (Bl3a) and (B13b) to simplify the expressions 
for the slope and the radius of curvature of the streamlines 


is as follows: 


(y+1) 72 o' Io re. s! 
a! o,=e uae (Oe Ae) {1 + £ os } ; 
A ij 20.. exp ( 2012 _ = 20 
i>l, 521 (C14) 


where it has been assumed throughout that the initial 


streamtube widths are uniform; that is, An! = lor 


1,j 
Any = 0 for all j. | 

An expression for a streamtube, which is far enough from 
the heat release region that it can be considered isentropic, 
is now derived. Carrying out the integration over y' in 
Eq. (Cl2a) and assuming that ni > 40' an expression for 


the slope of the bounding streamtube is obtained: 


dn: on (yt+1)Ijav27 oF s!? 
g — aon. e * « oo rs 
rJ : oI 
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Integrating Eq. (C15) with respect to ds; gives an explicit 


equation for the bounding streamtube perturbation from its 














5 initial position nj} ., as: 
1,3 
‘ see] 
pe es ee aes 
i,j* "1,5 ° 2 ce. ANG Q41 
(y#1) I a0"? Si | | 
= 10 (l + erf=——) . (C16) 
(2) 20° 
1 
ig 
f 
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APPENDIX D 
INTEGRAL EQUATION DERIVATION FOR TREATING 


SMALL DISTURBANCE TRANSONIC FLOW WITH SHOCKS 
FOR AN INFINITE WALL WITH GAUSSIAN SLOPE 


1. GENERAL DERIVATION OF INTEGRAL EQUATION 

7 In this appendix, an integral equation is derived for 
the small disturbance transonic flow with shocks on an 
infinite wall with Gaussian slope. The small perturbation, 
steady, two-dimensional transonic equation without heat | 


addition can be written as: 


ae uke 25 2 7 
(LM) drys Hoyrys = CLM ML (YtL) O51) Ogi ge Hoytys = 0 
(D1) 


where the local Mach number is approximated by the freestream 
Mach number plus a perturbation quantity (M7 (y+) O04) « 
Equation (Dl) is valid only in regions where the necessary 
derivatives exist and are continuous. Therefore, the integral 
equation derivation includes an equation for the transition 
through a shock; the equation is the classical relation for 


the shock polar [Ref. 34], Eq. (D2). 


(D2) 
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where subscripts a and b refer to conditions ahead of and. 
behind the shock. If a small perturbation analysis is 
carried out on Eq. (D2) similar to that performed in 
Shestaing Be: (D1), the following relation is fourd for 
the perturbation velocity components across the shock wave: — 
© (lem?) (ap = BE) + (HE BY? = MLA evan) (BE + OY) (HY - GY) 7/2 

| : (D3) 


Equation (D3) corresponds to the shock polar curve for 

shock waves of small strength that are inclined at any angle 

between that of normal shock waves and that of the Mach lines. 
The apolteabie boundary conditions for the integral 

equation analysis are that the perturbation velocity must 


vanish far upstream of the disturbance, Eq. (D4): 


= 0 (D4) 








and that the flow must be tangent to the wall surface, 


Eq. (D5): 
a) saat 2 he i 
Byr ar gt ee arr) ee 
17J ] 


where 92Z'/8x' is the slope of the wall surface in the x' (s') 
direction which can be approximated by Eq. (C15); and T is an 


equivalent thickness ratio for the surface defined as the 
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. 


maximum wall growth t (cm) in the y-direction divided by a 


characteristic length /27 o for M = 1.0. 


— (y-1) I jane? 
Tao | | YP ote | 





=e 


Also, it is necessary to prescribe that the direct influence 
of a disturbance in the supersonic region proceeds only--in 
the downstream direction. | 

Sprieter and Alksne [Ref. 10] derive in detail the 
integral equation by a Green's function analysis for transonic 
flow around bodies having subsonic freestream velocities, and 
a local shock discontinuity. The nonlinear term in the 
differential equation was maintained. This derivation can by 
extended to the streamtube configuration of the boundary 
condition for the internal region. For the sake of brevity, 
only the important equations and pertinent aseuasetens will 
be given in the following analysis. 

Since the object of the analysis is to determine the 
pressure perturbation and/or velocity perturbation on the 
bounding streamtube, an equivalent equation for the velocity 
perturbation can be readily obtained by differentiating 
Eq. (D1) with respect to x'; it is: 


2 wy 2 
= 2) x1 ~ 4 = 2 d u' . 
(1 Meo )ULiys + yiy' ~ M (y+) mee (—- ). (D7) 
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For the Green's function analysis, it is advantageous to 


normalize Eq. (D7) with the following linear transformation: 





x Ext > yi By’ 7 
m2 (y+1) a! owt yt a | 
u = ee 7 VS c} (D8) 
Bo - B 
where 8 = A, =u? 
In this'‘way, Eq. (D7) reduces to the following: 
= 2- 2— —2 | 

a + 23 =y-y = -_ (S) (D9) 
ax oy 





where u can further be defined by the approximation to the 


local Mach number, Eq. (D1), namely m2 = M2 (1+ (y+1) a), as: 


M7(y+1)a" mM? - mM? 
ee ee 1 --M_” 


From Eq. (D10), it is clear that for a subsonic freestream 


Mach number: u < 1 when the local flow is subsonic, U = 1 





when it is sonic, and u > 1 when it is supersonic. 
The Green's function analysis is a Neuman problem for 
finding a function u in a bounded region R, bounded by a 
é finite number of closed smooth contours C, such that Eq. (D9) 
is satisfied in R, u and its first partial derivative are 
continuous in R (shocks are excluded by judicious selection 
of the bounding contour) and the normal derivative 9u/0y 


. is prescribed on the boundary of R, Fig. 1D. 
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FIGURE 1D. Integration Region for Green's Function Analysis 
on Flows Past Bounding Streamtubes at 
_ Transonic Speeds 


Green's theorem states that: 


f (ave-n-cvi-n) dc = ff av2c -cv7a) aR (D11) 


( 
Cc R 


where the directional derivative on the left side are 


s taken along the normal, drawn inward, to the curve Ce 
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€ 


The Green's function G is constructed so that it satisfies 


the Laplace equation yc = § (E-x,n-Yt¥,) with 9G/an = 0 on Cy: 


ai _ il " 
ErnrXrY-Y_) = oe oS 


i 
|r 
war) 
3 
2 


(D12) 


where 





(R-E)* + (¥-¥,-7)° 
| (x-F)? + (F+¥,+H) ° 


and = and yj are the variables of integration while * and y 





are the coordinates of a. point P. With the Green's function 


of Eq. (D12), Eq. (D11) becomes: 


i702, Eh ox 
W(%,¥) - ¢ 2 (3) ana 
_ yee 
= f one) ae- f (oS -4)  af 
+f Ra ae f oe B- F8h ac wi 
b on an xX, Co on on 


where the respective integration regions and contours are 
designated in Fig. 1D. Integrating the remaining surface 
integral by parts twice with respect to =, applying the 
Green's function boundary condition on C, (n=0 or y=y,), 


allowing the radius r of the curve C, go to infinity, and 
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decomposing the line integral over the shock wave into 


components parallel to the axes of the coordinate system 


yields Eq. (D14). 
: wipe en De = = Ls . =2 °C tan 
u(x,y) - 5 U(x,Y) = Ua + F j f gE 
-f{ is 3 (G-5 a) - G@-5 0) | 
i | ob ja 
= c 2a -507) - a- Fa’) 4 
ae d6 1b 
* + [ 5 2s | = c mag 26 sos (nen) } dn 
: on ania an on jb cos (n,&) 
: | (D14) 





Lw 
wall shape is specified. 


where U. represents the linear wall velocity u when the 


a= f ¢S-73 at (D15) 
«00 an on n=0 


Sprieter and Alksne [Ref.10] in their derivation delineate 
- -. the advantages for integrating the surface integral over R 
by parts as was done in arriving at Eq. (Dl4). The advantages 
are important enough that they are reiterated here. First, 


the double integral of Eq. (D13) shows a very strong influence 


80 











% 


of the velocities in the region immediately surrounding the 
point P since they are multiplied by 1ln(R). This influence 


isi largely nullified in the double integral of Eq. (D14) 


because part of the region has a negative influence and part 


has a positive influence. The predominant influence in the | 


latter case is furnished by the term su standing outside 
the integral. Next, the contribution of distant regions 
is diminished in importance in the double integral of 


Eq. (D14) since their influence varies inversely with the 


‘square of the distance, rather than with the logarithm in 


Eq. (Dll). A further advantage is that the value of the 
double integral of Eq. (D14) is continuous through a shock 
wave rather than discontinuous as is the case with Eq. (D13). 
Lastly, a point of great importance in the approximate 
solution arises from the fact that the integration by parts 
provides extra terms (those containing = 0) in the 
integrals along the shock surface S which combine with those 
already present in such a manner that the contribution of 
these integrals becomes very small when the shock wave 
approaches a normal wave, as is usually the case at high 
subsonic speeds. 

Therefore, u must satisfy Eq. (D14) where U,, is given 
by Eq. (D115) and the shock relationship Eq. (D16), which is 


the normalized form of Eq. (D3). 


2 


(G,-O)? + H-¥)? = 5 G@ +0) @,-F) (D16) 
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In the case of a specified symmetric wall shape, Eq. (D14) 
can be considerably simplified as it permits the introduction 


of the boundary conditions which are specified at the outset 





“se 
of the problem, Eqs. (D17a) and (D17b), into the integral 
° over the wall surface. 
u(x,0) = Urw = KNOWN (D17a) 
and 
ss ae ns 2— 2 ~~ 
Oy a 88g 28) we es ’ (D17b) 
ott 3 OG ag ag 20 | 
where Z represents the reduced ordinate of the surface, and 
7 t the reduced thickness ratio given by: 
M7 (y+1)t 
: —— (D18) 
B 
Substituting the relationships of Eqs. (D17a) and (D17b) into 
Eq. (D115) completely determines at the outset of the analysis 
the linear wall velocity distribution: 
es ae : > x 2= 7 
Up (%,0) = _lim ( = f (1n (=) <2) dé 
Y"Yo6 aed 2 dt 7=0 
co — 
1 az d 
-bs & (019) 
-~ d& (x-€) : 
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2. SIMPLIFICATION AND APPROXIMATE SOLUTION OF THE 
INTEGRAL EQUATION FOR TRANSONIC FLOW 


To obtain a solution to Eq. (D14), some further 


assumptions must be made. They are: (a) all shock waves 


‘lie in a plane transverse to the beam, and (b) the shock 


waves are normal (i.e. normal to the local flow direction). 
These two assumptions reduce Eq. (D16) to an expression, 
Eq. (D20), which permits an advantageous introduction into 


Eq. (D14) to eliminate the integral over the shock S, Eq. (D21). 


(G-5 0), = (@-5 U), 
(D20) 
2 Pt Sey ees dae eS Se 
—— (u-s ue). = — (u-F u") 
coo 2 =. 2 rz | 
WE) = EGY) +7 PY + oe Sf S Sy inte) anak 
-20) 0& 2 (D21) 


Equation (D21) will be used as the basis for obtaining 
the velocity distribution on the bounding streamtube. 
Approximate solutions of this equation are obtained numerically 
by starting with an assumed velocity distribution and 
iterating until convergence is obtained. If mixed flow 
exists, the initial velocity distribution must include a 
proper discontinuity complying with the shock relations of 


Eq. (D20). 
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Due to the double integral in Eq. (D21), an iteration 
process would be cumbersome unless it could be reduced to 
a single integral by introducing suitable approximations. 
Oawatitsch [Refs. 35 and 36] has shown that approximate 
knowledge of the velocity distribution in the vicinity of 
the surface can be expressed in terms of the local coor= 
inate y, the bedinates of the wall shape Z(x), and the 
desired but unknown velocity distribution G,, (X, 0) on the 
surface. This permits the reduction of the double integral 
to a single integral. 

Oswatitsch has considered such an acoroniaate relation 
in which the velocity U(x,y) starts from the value u(x, 0) 
at the surface with an initial rate of change given by the 
irrotationality condition, Eq. (D22), and vanishing at large 


distances as 1/7. 








sup = (D22) 
ay Se ox w 

The resulting relation takes the following form: 
ee a, (x, 0) | 
u(x,y) = ———____—__, (D23) 


(L + (¥-¥,)/b) 


where b is a function of x satisfying the irrotationality 
condition at y = Yo° Solving for b using the boundary 


conditions of Eqs. (D17a) and (D17b) yields: 
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= b(x) = - —~—,- . (D24) 
Pag : d 


Substituting Eq. (D23) into the double Phesavad of Eq. (D21), 
integrating with respect to nh and setting y = Yo results in 

a simplified integral equation for calculating u,,(X, 0). 

For computational purposes, the integral equation is expressed 


ass 





- — _- Land paw ae: 
WW, (0) = Wry (Xe 9) + Zu, (x, 0) 51 (D25) 
where 
co a (E,0) F- ame 
=O : 
and the function E is: 
(I +X") 
(1 - 10x? + 5x4) in| x| - (14x) (25-71x7-x4-x9) 


Although I is a function of u, and is-unknown, it is informative 





to rewrite Eq. (D25) by solving for a, in terms of I and ran 
thus: 


) . ual ce (2u,,71) =1 s/I - L , (D27) 
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where 


L 


2uL - il 


a known function since Mie 


specified. It should be noted from Eq. (D27) that the 


is known once the wall shape is 


discriminant must always be positive in order to obtain real 
values for Uo thus I >L. Furthermore, the choice of the 
plus or minus sign determines whether the local velocities 
are subsonic or supersonic. A change in sign at a point 
where the radical is zero corresponds to a smooth transition 
through the sonic speed. A change in sign at a point where 
the radical is not zero corresponds to a discontinuous jump 
in velocity, namely, a normal shock from supersonic to sub- 
sonic flow when progressing in the flow direction. A dis- 
continuity in the reverse direction is inadmissable since 

it corresponds to an expansion shock, an impossible phenomenon 


which violates the second law of thermodynamics. 


3. COMPUTATIONAL TECHNIQUE 

In this section, a brief description of the computational 
technique used to solve Eq. (D25) will be presented. Sprieter 
and Alksne [Ref. 10] describes in great detail how a judicious 
selection of the initial velocity distribution, depending 
upon the surface shape and flow configuration (i.e. entirely 
subsonic flow or mixed flow with possible shocks), is essen- | 


tial in providing a rapid convergence to the solution. It 
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will be assumed for the following discussion, therefore, 
that a judiciously selected velocity distribution is known. 
The method of solution is as follows: 

a. From the known surface shape in normalized coor- 
dinates Z, calculate its first two derivatives with respect 
to the normalized flow coordinate x; that is, dZ/dx and 
a*Z/ax-. Since the normalization depends upon the freestream 


Mach number, which is unknown, the value of T, Eq. (D18), 


will be a parameter in the iteration technique. 


b. Calculate the linear wall velocity Wy py (0) from 
Eq. (D19). As T does not depend on €, it can be taken 
outside of the integral of Eq. (D19) giving an expression 


directly proportional to the parameter T. 
Uy py (%20) = T Wh (x, 0) (D28) 


where Wry 6%» 0) is a function of only the position x and is 
fixed once the surface shape is specified. The quantity 
L follows impadiately from Eq. (D27). 

c. Assume an initial value of T consistent with the 
velocity distribution and calculate b(x) from Eq. (D24). 


d. Since b(x) is now known and an initial wall velocity 


(u,(x,0) has been selected, Eq. (D26) can be integrated 


numerically to get I. 
-@. The I function calculated for the initial value of 
Tt is graphed. T is then varied in the L function until. 


the two curves are tangent at, at least, one point (i.e. 
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I-L = 0). As was stated earlier, I > L at all points x 
nigneceie surface for the velocity to be real. 
| f. Using the L function that satisfies the tangency 
requirement, a new approximation to the solution is obtained 
from Eq. (D27), where the sign chosen for the radical is 
important. | | 
g. The iteration procedure continues with the new 
a, (x, 0) and T until both quantities converge to a solution. 
h. Once a solution is obtained, the freestream Mach 
number for a given flow condition is calculated from 
Eq. (D18) and the local Mach number distribution from 


Eq. (D10). 
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APPENDIX E 


HODOGRAPH TECHNIQUES FOR TRANSONIC FLOW 


The hodograph plane represents a coordinate transfor- 
nection in which the velocity components are the independent 
.. vwariables. Spatial coordinates x and y are replaced by | 
velocity components t and v; that is, the physical esordi= 
nates are represented by x(i,¥) and y(i,¥%) [Ref. 26]. This 
. coordinate transformation is convenient because it enables 
the simple presentation of data or solutions and, more 
importantly, it linearizes the two-dimensional planar tran- 
sonic flow equation that is nonlinear in the physical plane. 
This equation is known as the Tricomi or Tricomi-Euler Equa- 
tion. Boundary conditions are difficult to satisfy [Ref. 15]. 

The small perturbation transonic equation without heat 


addition can be written as: 


42,90" | ae! _ 2 rt 
(1 M,.) 3xT + Byr = Mo(ytl)u xr (El) 


hoa | wy t 
gor - ger = 0 (E2) 


Using the hodograph transformation, Eqs. (El) and (E2) | 


become: 
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(1-2) OY) 4 2X" = wt (yeayur SO (E3) 
av oat ae! 


and 


ox Mag | (E4) 


ei (E3) and (E4) are the linear transonic hodograph 
equations. | 

As is mentioned in the introduction, the linearized 
small perturbation equations for supersonic flow are hyper- 
bolic and those for subsonic flow are elliptic. The 
transonic hodograph equations change from elliptic to _ 


hyperbolic when the coefficient of dy'/av' changes sign: 
2 2 i. 8 
(1-M) - Mo(yt1)u' = 0 (E5) 


The critical velocity, where the change of sign occurs, 


can be directly solved from Eq. (E5). 


(ew?) 


u'* egy 3 (E6) 
(y+1) My . 


Therefore, for ut < u'*, Eq. (E3) is elliptic, and when 

u' > u'*, it is hyperbolic. The features of transonic flow 
are captured by the transonic hodograph equations. Since 
the flow occuring in the external region and on the bounding 


streamtube is considered to be isentropic two-dimensional 
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and planar, the hodograph method gives a conceptual idea of 
what the resulting flow configuration should be over an 
infinite wedge of prescribed shape. 

Guderley and Yoshihara [Ref. 15] use the flow angle 6 
and n, defined in Eq. (E7), as the independent variables. 


1/3 ¢ (yewh) Jw) (E7) 


n= (y¥+1) 
where w= fu? + v? - With this selection of independent 
variables, the streamfunction in the hodograph plane is: 


Van 7 1 Vee = 0 (E8) 


nn 


and the velocity potential function is: 


ae — TN dg9"= 0 (E9) 
The formulation is Soe a flow with constant stagnation 
temperature and, hence, can be applied to the solution for 
the bounding streamtube. For flow over a wedge, Guderley 
and Yoshihara found a steady solution for M, precisely unity. 
Fig. lEa shows the flow over the bounding streamtube 
in the physical plane. Points A, B, C, D, and E are shown 
mapped in the hodograph plane, Fig. 1Eb. A limiting Mach 
line is enue in the physical plane. Any characteristic din 


the region bounded by the sonic line and the limiting Mach 
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(a) Flow in Physical Plane (Not to scale) 
\ Supersonic Flow 
Subsonic Floy \7 | 
C +Characteristics 
= | (b) Flow in the Hodograph Plane 
FIGURE 1E. Sonic Flow over Bounding Streamtubes 
i | (Physical and Hodograph Plane 
Representations) 
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line can interact with the sonic line. Any characteristic 
downstream of the limiting Mach line does not interact with 
* the sonic line, and, hence, the flow is not affected by the 


subsonic flow upstream of the sonic line. 
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APPENDIX F 


DERIVATION OF THE SLUING RATES NECESSARY TO 
PRECLUDE SIGNIFICANT PHASE DISTORTIONS IN 
A LASER BEAM AT MACH 1.0 


When the phase difference between two BGiage in a laser 
beam is gyéater than a quarter of a wavelength, the deviation 
of the beam from its original direction is sufficient to 
‘cause considerable distortion. The range of laser sluing 
rates necessary to prevent such phase distortions and, hence, 
to preclude defocusing of the beam at Mach 1.0 is derived. 

Consider a laser beam that is slued as shown in Fig. l. 
The Mach number of the air relative to the beam varies along 


the beam as: 
= eee (Fl) S 


All Mach number regimes occur along the beam. If L repre- 
sents the radial distance between two stations on the beam 
where the Mach numbers are M, and Mo, it can be determined 
from Eq. (Fl) that L is equal to: 


(M,-M,)a 
a. (F2) 


The number of waves (N) of laser radiation in the region 


of length L is given by: 
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“> | _ Ns += a (F3) 


where A is the wavelength of the laser beam and is equal 
to nv/c. v (radians sec +) is the laser frequency equal 


13 Hertz for a CO, laser, and c (m 


to approximately 2 x 10 
sec"). is the speed of light. Therefore, the change in the — 
number of waves in the beam length ts, due to a change in 
the index of refraction, is obtained by differentiating 


Eq. (F3). 
~ | | an = vidR | (F4) 
The change in the refractive index is related to the 


density perturbations in the beam, a known quantity, through 


the Gladstone-Dale Law [Ref. 37], a simplified version of 





the Lorentz-Lorenz Law for gases with n = 1, Eq. (F5): 
n=l+t+ko, (F5): 


where k is the Gladstone-Dale constant which is approximately 
equal to 2.25 x 107? mikgm for air with wavelengths beyond 
the near infrared (\ > lu). Differentiating Eq. (F5) and 
substituting it and Eq. (F2) into Eq. (F4), the desired 
expression for the change in the number of waves is obtained 


as: 


k(M,~-M,) avo, 
: | pee oe (F6) 
ae : Qc Poo 
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It is of interest to determine the range of sluing 
rates for which the change in the number of waves is less 
than a quarter of a wavelength at Mach 1.0 from Eq. (F7): 

4k (M 27M,) ave,, 


i ne“ 


c Max — Puin) : Ev) 
where Omax 25 the maximum positive density perturbation 

and Blin is the minimum (maximum negative) density perturba- 
tion in the beam at Mach 1.0, points A and B in Fig. 1l, 
respectively. As an example, using M, = 1.001, Ms = 0.999, 
the flow properties in Table I, and the maximum difference 
in the density perturbations at Mach 1.0 from Fig. 11 of 
1.5 x 1074, the sluing rate must be greater than 7.5 x 1073 


radians/sec. 
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APPENDIX G 


SUBSONIC AND SUPERSONIC THERMAL BLOOMING 
AND COMPUTER PROGRAM (BLOOM) 


Tsien and Beilock [Ref. 38] have obtained solutions for 
the linearized equations of motion with heat addition for 
subsonic and supersonic flow. The solutions are for a line 
heat source q in a uniform, two-dimensional planar flow 
which can be considered to be a perturbation on an initially 
one-dimensional flow. 

In 1971 and 1972 studies were conducted at NPS by Fuhs 
[Ref. 39] on the applicatin of external heat to the reduction 
of drag. These studies required thorough analysis of heat 
~ addition in both subsonic and sapesenhideriowas The analysis 
of the flow with heat addition found two other applications. 
One application.was external burning assisted projectile. This 
application was studied by CDR. W. J. H. Smithey [Ref. 40]. 
The other application was recognized by Fuhs [Ref. 4]; this 
involves iene density inhomogeneity arising from quantum 
inefficiency. Professors Biblarz and Fuhs [Ref. 5] developed 
the capability and extended the analysis to predict density 
inhomogeneity in a large gas dynamic laser; see the paper 
by Fuhs, Biblarz, Cawthra and Campbell [Ref. 6] for details 
of the experiment. Subsequent to the application of linear- 
ized solutions to flow internal to a laser, it became apparent 
that the analysis could predict flow properties in thermal 


blooming. In a series of papers and presentations, Fuhs, 
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Biblarz, Burden and Carey {Refs. 7, 8 and 41-43] developed 


the appropriate analysis for linearized thermal blooming. 


1. SUBSONIC FLOW WITH HEAT ADDITION 


The perturbations to the flow [Ref. 38] are given by: 


2 
it 


__y-1)g x | | | ~(g1) 
2TYp.,Pu,, (x2 + B2y2) 


~ (¥-1) Bq | 
ee | - (G2) 
a. (Y-1) Mg = 
Pp ™ ~ ZialBp,, (x? + B2y2) + B2y (G3) _ 
(y-1) Mog 


(y-1) 5 ( 
eo ~ocps y) I(x) (G4) 
(x? +B? y) Mco 


D 
-~ 
I 


where the symbols are defined in the list of symbols and 
6 = Vv71-mM* for subsonic Mach numbers. The second term in 
Eq. (G4) is the wake of the line heat source q. The 
temperature perturbation can be obtained from the equation 
of state. 

Biblarz and Fuhs [Ref. 6] in applying Eqs. (Gl) through 
(G4) to laser internal aerodynamics have presented the 
integral equations used to solve the unbounded heat release 
region problem (i.e. subsonic thermal. bi doming) and are 
listed below for completeness. Since the heat source is 
volume distributed, instead of a single infinitesimal rod, 
an integral representation was used to account for all 
contributions to the flow properties over the wOLuiie: The 


geometry involved in the calculation is shown in Fig. 1G. 
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FIGURE 1G. Geometry of Flow 
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The heat quantity h(x,y) is related to the laser output 
intensity through the kinetic lag that may be present. To 
develop an equation incorporating the lag between laser | 
energy and heating, it is necessary to examine the kinetics 
of relaxation. An exponential function describes the 
relaxation with a relaxation time t[Refs. 44 and 45]. The 


heat release function is given by: 


x ' eS 
a f £(x",y) exp[- 2-22) ax" (49) 


h(x,y) = 1 
Y o co fo <] 


where 0 e Too f(x,y) is the laser energy per time per 
volume. The volume can be interpreted as an area times the 
span of the heat addition zone across the planar flow. 

If significant variations of beam intensity in the radial 
(beam) direction exist, the two-dimensional heat addition 
function can be approximated from the three-dimensional 
h(x,y,r) by an average: 

| r 
hizy) wok: J h(x,y,r) ar 7 (G10) 
i | 7 


OF ee 


where Ar = 1r,,,-r, .- The beam is divided into segments Ar 
in length and Eqs. (G5) through (G8) are solved for each | 
segment. 

If the heat release h(x,y) is constant, the equations 
can be integrated in closed form. The easureine dateqvere 


are: 
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. | | u' = Sryphay M(xry7&, mMpB) (G11) 
; . ee (y-1)Ioo | 4 | 
id Vv" = Bryp,6u,, N(x,y7&,n,B) (ey) 
a | (y-L)M,Toa | 
p == Ta Bp, | NCTE nl (G13). 
| (y-1) M,, t,o 
p'=- Ky 2 M(x,y?&,n,8) 
Bp, 
(y-1)I,0 | 
——— ae I(x-&) dé&dn (G14). 
Ay UpP. 
m - where 
2 
Mix,yr&.n/8) = | | 25% an al 
7 | g? (y-n) 
id 149 bo 
| a (re) tan | 8 (y-n) 
(x-§) | n E 
. 1 


| 2 2 2 
N(x,y?&,n-8) = ap in es is (y=n) | 
28 B . 
No +52 
_ (x-&) (y-n) -1l]| (x-&) 
g* a 8 van eS | Ny aT 





The first term in Eqs. (G8) and (G14) correspond to the 
density perturbations in regions outside the wake and the 
last term represents the contribution of the wake and can 
be interpreted as the consequence of the heat convected 


‘ along the streamlines by the flow. 
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2. SUPERSONIC FLOW WITH HEAT ADDITION 
The perturbations to the flow are also formulated by 


v  ™gein and Beilock [Ref. 38] as: 


; | er. tye) a | | ve 
= i! = eeu, 8 By) | | (G15) 


(y-1) | 
Pee ox -~ By) | (G16) 


y' 


(y-1) M4 | 
3a tp. 6(x - By) (G17) 


S(x - By) - HOUS siy) r(x) (G18) 


Bp. M,, 200 Poo 


(y-1) Mg 
S 





oO 
- 
Nt 


2a,, 
‘ | Fuhs [Refs. 4 and 39] has modified these equations into_ 


an integral form as: 








s 
~ (Y-1) ; 
uf =- : h(é,n) sin u ds (G19) 
2Y P,P Us, : 
o = gli) f° h(E,n) sin u ds (G20) 
2° Poo 0 ‘ 
(Y-1)M, °° 
Bt = h(&,n) sin wu ds (G21 
P 28, PP oo : ’ | ) 
ee 8 ot Ot | 
pt = s— f h(E,n) sin u ds 
2a, BP, 0 
A af i AYES) oh 2 
Os f f n(én) S(y-n) I(x-&) déan (G22) 
: Mia Po §& 1 3 
- where the symbology is the same as that used in the subsonic 





equations except that 86 = Veet for supersonic Mach numbers. 
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3. COMPUTER PROGRAM i: 
The FORTRAN computer program BLOOM, included at the end 
of this appendix, was written to compute and/or plot the flow 


~ 


eee acai (a', v', 6' and 6') in the subsonic and super- 


sonic’ ‘low with heat addition. Using a modified version of the 





computer program DIDER [Ref. 5] as a subroutine, the subsonic 
perturbation equations, Eqs. (G5) through (G9), and the | 
supersonic perturbation equations, Eqs. (G19) through (G22), 
are numerically integrated to obtain the linear thermal 
blooming solutions. The program calculates the two- 
dimensional field at any cross section of the laser beam 
where the energy release distribution is known. Finite 
kinetics are incorporated as described by Eq. (G9). 

The program (BLOOM) is versatile in that it can calculate 
the flow quantities for several subsonic and/or supersonic 
freestream Mach numbers, heat release distributions, relaxation 
rates, and freestream flow properties in a single computer run. 
Since the length of time for each individual case depends 
primarily on the size of the flow ere and on the mesh that 
is superimposed upon it, time requirements must he determined 
on an individual basis. For the results presented 1 Figs. 4 
through 7 and Figs. 12 and 13 from the program BLOOM, the | 
two-dimensional region of interest was 60 cm. square centered 
at the origin of the Gaussian heat release distribution and 
was descretized into a mesh of cells l Sia square for 
subsonic flow (M, = 0.999) and into rectangular cells whose 


size is determined by the Mach angle for supersonic flow 
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(M. = 1.001). Supersonically, the ratio of cell height Ay 
66 cell length re is equal to tan yu. The subsonic results 
take approximately 25 minutes and 140 K of memory (without 
plotting subroutine CONTUR) while the supersonic results 
take approximately 3 minutes and 285 K of memory (without 
plotting routine CONTUR) to compile and execute using an 
IBM 360 computer. | 

The program can also calculate constant heat release 
distributions for single and multi-pass laser beam configura- 
tions with or without wall reflections and kinetics. In its 
present form, the program can not calculate supercritical 
(mixed) flows across the beam. 

A main feature of the supersonic portion of the program 
is that it sets up the characteristic lines belonging to each 
individual cell (i.e. their diagonals) and adds up the 
contributions for each flow property from the characteristics 
and, in the case of the density perturbation, from the wake 
at each cell. The calculations along the characteristics 
are facilitated by the choice of cell shape and size as 


described above. 
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